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Abstract 

We study the dynamical behavior of the Sherrington Kirkpatrick 
model. Thanks to the APE supercomputer we are able to analyze 
large lattice volumes, and to investigate the low T region. We present 
a new and precise determination of the remnant magnetization and of 
its time decay exponent, of the energy time decay exponent, and we 
discuss aging phenomena in the model. We exclude validity of naive 
aging, and propose different options that fit the numerical data. 



The study of the dynamical behavior of disordered (and complex) systems 
is receiving a large amount of attention (see for example |[ |3|, |j] and refer- 
ences therein). On one side this is because the study of dynamical behaviors 
can shed light on a large realm of new and interesting phenomena (aging 
is one of them). On the other side a reliable description of the amorphous 
glassy state will surely include a crucial, non-trivial dynamical behavior. 

In many cases a full analytic computation is not yet feasible and many 
analytical results are based on conjectural grounds. Numerical simulations 
are here very useful in order to support the theoretical conjectures. Un- 
fortunately a critical issue in the dynamics is the dependence of the times 
needed to reach equilibrium on the lattice size. This issue can be clarified 
only when both the size of the system is very large and the observation times 
are very long. This problem may be bypassed by considering the behavior of 
the infinite (practically very large) system as function of the (Monte Carlo) 
time. 

This approach requires the study of very large systems. Unfortunately in 
the models that have been studied better, for which a full analytic solution for 
the static behavior is available, the time needed for executing one dynamical 
step increases severely with the lattice size. It grows as the lattice size N to 
the second power, N 2 for the Sherrington-Kirkpatrick model (SK) [||], and as 
N p with p > 3 for the p-spin models. Faster simulations can be done for the 
diluted models, but in this case we do not known the static solution exactly. 
Moreover, given our lack of analytic control of the dynamic behavior, we 
cannot be completely sure that there are not some subtle differences among 
diluted and non diluted models. 

Long range models may be classified into two different categories: 

• In the first class of models (taken in the infinite volume limit) intensive 
quantities like the energy E(t) or the magnetization m(t) evolve toward 
their equilibrium values in the limit where t — > oo. 

• In the second category dynamical intensive quantities do not tend to 
their equilibrium values, and truly metastable states are present. 

Such difference in the dynamical behavior is strong. It is believed that 
the models with a continuous replica symmetry breaking belong to the first 
category while models with a one step replica symmetry breaking belong to 
the second category. 
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The present investigation addresses this question in the case of the SK 
model 0. The issue is very sensitive because some numerical investigations 
have suggested that the asymptotic properties of the SK model are different 
of the equilibrium one and in particular that the remnant magnetization 
is non zero also at infinite time |7| (that would also imply that the weak 
ergodicity breaking scenario cannot hold). One of the results of this note 
is to show that the value of the remnant magnetization is compatible with 
being zero in the infinite volume limit and that the SK model does belong 
to the first category we have described. 

In the following we will show and discuss some of our numerical results. 
The main features that make these results relevant are that we have been 
able to work on large lattices and for low value of the temperature T. Our use 
of the APE- 100 supercomputer |J has been crucial to allow such large scale 
simulations. The fact that our program is truly parallelized (by dividing the 
spatial lattice among the different processors) makes possible to study very 
large lattices. Anyway the dynamics is still a sequential Metropolis one (we 
change one spin of a given system at the time). We also notice that since 
we are dealing with the infinite range mean field model obtaining, as in the 
code we usedQ, an effective parallelization on a mesh with fixed connectivity 
is far from trivial. 

In this note we will focus on the four main points we have been able to 
analyze: 

1. We determine with good precision the time decay exponent of the mag- 
netization. 

2. We analyze the behavior of the remnant magnetization as a function 
of N. 

3. We analyze the energy time decay exponent. 

4. We discuss aging through the magnetization-magnetization time de- 
pendent correlation function. 

First a few details about our simulations. We study the usual Sherrington- 
Kirkpatrick model 0, 0], with quenched couplings J chosen from a Gaussian 

lr The code we have used for the numerical simulation is due to P. Paolucci and D. 
Rossetti, unpublished. 
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distribution. We have studied lattices of size N = 2088, 4032, 8064 and 
18432, and T = 0.25, 0.50 and 0.75 with a local Metropolis dynamics starting 
from a fully magnetized state (all spins set to +1 state). We study from 100 to 
200 realizations of the quenched couplings for each N value but for the largest 
one, where we have order 15 samples. In each realization of the quenched 
disorder we have followed six copies of the system evolving independently 
(mainly to get a better computational efficiency). Runs have been 10000 
lattice sweeps long. 

To fit our data we have used the Minuit library from Cernlib, and the 
jack-knife approach to compute statistical errors. The results we present here 
have been obtained by using large time windows for the fits. 

The dynamical behavior of physical quantities in systems with relevant 
quenched disorder is usually expected to converge to the equilibrium values 
with power laws. The magnetization per spin on a lattice of size N is expected 
to behave as 

mW^m^+mfV" 1 , (1) 

where the superscript (JV) (that we will omit in the following when we can 
do so without creating ambiguities) indicates that the parameters depend on 
the lattice size. 

The previous formula needs a few comments. In any finite system the 
residual magnetization must eventually go to zero at infinite time. However 
we can distinguish at least two time scales. Let us consider for example the 
case of an unfrustrated system with two equal free energy states differing by a 
global spin reversal. At first the system will fall in one of the two equilibrium 
states, and on a much longer time scale the system will start to oscillate 
among the two different equilibrium states. Only in this second region of the 
time evolution the magnetization will become zero. If the ground states of 
the system are in some sense random, we may expect that from the central 
limit theorem the residual magnetization at the end of the first phase of the 
dynamics will be proportional to iV~2. A similar phenomenon for large, but 
not too large times is expected to take place here. The effective remnant 
magnetization is likely to be N -1 ^. While in an unfrustrated system uj — |, 
the presence of a large number of equilibrium states suggests a smaller value 
of uj and our data are compatible with the value uj = j found close to T c for 
the SK model. The possibility uj = 0, i.e. of a non decay of the remnant 
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magnetization with the size, is a priori possible, but it is not the one preferred 
by our data which are well compatible with an approach to equilibrium for 
all the quantities. 

In figure (|]) we plot the exponent S as a function of T for the different iV 
values. Our data are not precise enough to detect a clear N dependence of 
5 (only for the smallest lattice, N = 2088 one can maybe read a systematic 
deviation). The straight line that goes through the figure is the line 6 = 
it fits very well the data, and indeed the fit is best on the larger lattice size, 
N = 18432. We consider this plot as good evidence that 

W ^ £ • (2) 

c 



The exponent 5 is linear in |r, and 

J- c 



S (T - T") = 1 (3) 

(for a discussion of the difference of this limit value and the limit S(T — > 
T c + ) see the second of [0). Our estimate is also compatible with the best 



estimate of [10], while the models defined on 3 graphs (that is expected to 



have the same critical behavior of the SK model) seems to prefer a value 

8(F) * € Q- 

We plot our best estimate for the remnant magnetization as a function 
of T for different iV values in figure (H). The N dependence of the data is 
very weak. When assuming a power decay one finds that the data are fully 
compatible with the power |, as discussed in the second of IJ. It is clear from 
our data that a constant behavior, with a non-zero residual magnetization, 
cannot be excluded. 

Our third result concerns the energy time decay exponent. In figure (y) 
we show it versus the inverse square root of the volume. At low T we get 
an iV-independent estimate, with a clear dependence over T (we get .4 at 
T = .25 and .6 at T = .5). At T = .75 we seem to have a strongly iV- 
dependent result. Even if in this case it is not easy to be sure, the behavior 
of the power exponent is again compatible with a linear dependence on T. 

Our last results concern aging. We have measured the two time spin-spin 
correlation function 



C(t w ,t w + t) = jjQ2(Ti(jt w )ai(t w + t)) , (4) 
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always starting from random initial conditions. Our simulations were per- 
formed at T = .5T C , with "lattice size" N = 4096 and 16 independent samples 
with 6 replicas each. 

The simplest scenario one can think about is naive aging, i.e. 

C(t w ,t w + t) = f(U , (5) 

in the region where both t and t w are large. 

Our results for C are displayed in figure @. Although the data are in 
some very rough agreement with a naive aging behavior there are strong 
systematic corrections which obviously modify the functional form of (|5j). 
Naive aging is not satisfied in the Sherrington-Kirkpatrick model. The type 
of violations of naive aging are similar to those observed in real experiments: 
the function C(t w ,t w + t) at fixed s = j- decreases (increases) when t w 
increases for s < 1 (for s > 1). 

In naive aging one assumes that the time scale of relaxation at t w scales 



as t w and this assumption is not correct in the approach of |T2] and [g, [J. 
One can take one of many different attitudes. 

The first possibility is that naive aging is not correct, and we decide to 
use some phenomenological corrections to it. In this different form of 

scaling may be expected. The most popular assumption is interrupted aging 



[O. It corresponds to assume that 



f l+M 

C(t w ,t w + t) = f(- ). (6) 

The introduction of the power correction [i may help to fit the data at large 
j- but it does not improve the situation at small and spoils the agreement 
of aging at j- — 1. It seems that this correction is not very useful here. 
Moreover there is no theoretical justification for interrupted aging in this 
context. 

Another possibility, natural in the contest of @, |3| is to assume that 

C(t w ,t w + t) = f( l ^±J^l h (7 ) 

We have plotted our data for C(t w , t w + t) versus ln(t + t w )/ ln(t w ) in figure 
(151). This change of the scaling law definitely improves the situation, so this 
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is a possible solution to the problem. It is interesting to note that in this 
version of the scaling form we have that 

lim C(t w ,Xt w ) (8) 

t w — >oo 

goes to a value independent on A, which can be therefore identified with q^A- 
Obviously in this case the function / will be discontinuous when its argument 
is equal to one. We can also say that the time for reaching a given value of 
C < qEA is given by 

t = tt C) (9) 

where the function 7 is related to the function /. For large values of t w the 
value of C that can be reached after the time t and after At do coincide. 
In this way the ultrametricity of the configurations is satisfied also in the 
dynamics (as it should be): for example it implies that if in a given time we 
can go at distance 1 — C from a given configurations, we can arrive always 
at the same distance if we double the time. This would be the scenario that 
is in better agreement with the existing theoretical computations 0. 

The last possibility is that naive aging is asymptotically correct, but 
there is a correction to it which vanishes as a negative power of time. The 
rational for this choice is that we know (as we shall see later) that some 
power corrections to scaling are necessarily present at small t. 

Here we explore if this third possibility is compatible with our numerical 
data. We note that in the region where t <C t w we should have 

C(t w ,t w + t)=q EA + At- a . (10) 

The exponent a is a function of T; it is equal to 0.5 at the critical temperature 
and it decreases to a smaller value when T decreases (it has been estimated 
to be equal to 0.36 at zero temperature As far as the value of a is not 

too large (and it is certainly much smaller than 0.36 for three dimensional 
spin glasses, i.e. close to 0.06) the corrections to the scaling law cannot be 
easily neglected. 

The simplest modification to the naive aging prediction which is compat- 
ible with the previous equation is: 

c(t w ,t w + t)=f (U+r a f 1 (U, (11) 
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where the two functions fo(s) and fi(s) are not divergent in the limit s — > 0. 
In principle also the exponent a could be a function of s, but for simplicity 
we assume that it is a constant. 

We have fitted the data using equation ([□]). The best value of a we 
find is 0.2, which is a factor two smaller that the theoretical predictions. 
The origin of this discrepancy is not clear (possible reasons are finite volume 
effects, a crossover behavior at s = 0, a strong s dependence of a). If we stick 
to our best value a = 0.2 we obtain the functions fo(s) and fi(s) displayed 
in figures @ and (0). In this way we estimate a value of qEA close to 0.75, 
which is slightly different from the theoretical value ~ 0.6. Also in this case 
the origins of this discrepancy are not clear. 

In order to exhibit the quality of our best fits we show in figure (|8p the 
quantity C(t w ,t w + t)— fi(r-) t~ a . Here, at the expense of having introduced 
an extra function, the data seem to collapse well on a single curve. We also 
plot the same quantity for a different value of a, to show that in this case 
the collapse is not as good. 

The main difference among the behavior of ( [TTj) and the one of (|7|) is that 
in the case of ( |i~l~D lim^-^oo C(t w , t w ■ (1 + s)) is a non-trivial function of s, 
while in the case of (0) it does not depend on s. The scaling (0) is indeed 
in agreement with a picture where the barriers for reaching a value q of the 
overlap strongly depend on q. 

Longer runs on larger lattices will be able to improve our understanding 
of the situation. We believe that understanding details of the aging pattern 
is important, and that these results are a first step toward this interesting 
goal. 
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;ure 5: In C{t w , t w + t) as a function of ln(t + t w )/ ln(t„,). T = 0.5T C . 
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Figure 6: Fitted fo as a function of j- for a = 0.2. T = 0.5T C . 
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Figure 7: Fitted fi as a function of j- for a = 0.2. T = 0.5T C . 
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